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Abstract 



We establish rigorous results about the Newtonian limit of general rela- 
tivity by applying to it the theory of different time scales for nonlinear partial 
differential equations as developed in [4,1,8]. Roughly speaking we obtain a 
priori estimates for solutions to Einstein equations, an intermediate, but fun- 
damental, step to show that given a Newtonian solution there exist continuous 
one-parameter families of solutions to the full Einstein's equations -the pa- 
rameter being the inverse of the speed of light- which for a finite amount 
of time are close to the Newtonian solution. These one-parameter families 
are chosen via an initialization procedure applied to the initial data for the 
general relativistic solutions. This procedure allows one to choose the initial 
data in such a way as to obtain a relativistic solution close to the Newtonian 
solution in any a priori given Sobolev norm. In some intuitive sense these 
relativistic solutions, by being close to the Newtonian one, have little extra 
radiation content (although, actually, this should be so only in the case of the 
characteristic initial data formulation along future directed light cones). 

Our results are local, in the sense that they do not include the treatment 
of asymptotic regions; global results are admittedly very important -in par- 
ticular they would say how differentiable the solutions are with respect to 
the parameter-, but their treatment would involve the handling of tools even 
more technical than the ones used here. On the other hand, this local theory 
is all what is needed for most problems of practical numerical computation. 
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I. INTRODUCTION 



As it has been suggested through the extended history of the treatment of the newtonian 
limit of relativity, slow-motion corrections should be naturally obtained by suitable approxi- 
mation schemes starting from the newtonian gravitation theory. These schemes assume, for 
their validity, that there is a sufficient number of sufficiently smooth one-parameter families 
of solutions to the full Einstein equations (from now on referred to alternatively as the rel- 
ativistic equations) , the parameter being the inverse of the speed of light, and smoothness 
of these solutions with respect to this parameter being valid even at its zero value. Their 
smoothness is required so that the different order corrections can be interpreted as Taylor 
series coefficients of these one-parameter families. Their existence in sufficient number is 
required so that there is at least one for each newtonian solution, that is, so that we can 
describe every physical situation that can occur. Our purpose is to show some rigorous re- 
sults about this matter, that is to show the existence of a sufficiently large number of these 
smooth one-parameter families of solutions. 

The general theory of gravity coupled with matter sources obeying their own symmetric 
hyperbolic equations of motion is a problem of at least two different time scales, one having 
as characteristic speed the speed of light, the other -much slower- the speed of sound of the 
matter. One applies newtonian and post-newtonian approximations when one is interested 
in solutions where things happening at the fastest time scales are small and unimportant for 
the bulk motion of the sources. We will call them slow-motion solutions. That is the case 
for instance of the solar system; as a first approximation one is not interested in the details 
of the gravitational radiation, nor is this radiation important for knowing the motion of the 
planets in that approximation; we claim we can describe it -locally- with a slow solution. 
These slow solutions are the ones forming the one-parameter families referred to above. 

In recent years a complete theory has been developed for treating nonlinear partial 
differential equations with different time scales [4,1,8], which answers our questions in that 
general setting; not only does this theory tell us about the existence of these slow-motion 
solutions, but it also gives us a recipe called initialization on how to choose initial data for 
them. Even more, that recipe also allows us to choose data whose evolution is arbitrarily 
close to the corresponding solution to the limiting system, that is the system of equations 
obtained when the largest speed is set to infinity. Our approach to the Newtonian limit will 
then be to set the relativistic equations in a form suited for the application of the different 
time-scales theory above mentioned, apply it and see what kind of initialization procedure 
(further constraints on the initial data) we obtain. Solutions of this extra constrained system 
are the slow-motion solutions we are seeking. 

The plan of the paper is as follows: 

In the second section we present an introduction to the theory of quasilinear symmetric 
hyperbolic equations with different time scales. We briefly discuss what the basic principles 
of that theory are, what the requirements on the system of equations are for it to work, and 
finally what the results are. 

In the third section we present a new formulation (essentially a redefinition of variables) 
of Einstein equations as a first order symmetric hyperbolic system for arbitrary -but given- 
lapse and shift tensors. Here we assume the matter fields and their interactions with gravity 
to be such that the whole set of equations (including the matter fields) is block-diagonal 
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symmetric hyperbolic. For instance, the equations for the matter fields may themselves be 
symmetric hyperbolic and depend only on the metric and its connection (and not on the 
curvature tensor). The example we always have in mind as a matter source is that of a 
perfect fluid. 

In the fourth section we make use of this new formulation, i.e. of the freedom of choosing 
arbitrarily the lapse and shift, to pick a particular gauge in which our system satisfies the 
requirements of the different time-scale theory. This forces us to introduce elliptic equations 
on each time slice for the lapse and shift, and therefore to treat now a mixed hyperbolic- 
elliptic system of nonlinear equations. 

In the fifth section we discuss the initialization procedure, that is the selection of initial 
data that give rise to time-regular solutions. 

In the sixth section we summarize the results obtained and discuss how they might be 
embedded in particular settings to yield actual theorems. We also give our expectation for 
the local problem and for the asymptotically flat (global) case. 

II. THE THEORY OF PDE WITH DIFFERENT TIME SCALES 

There follows a short review of the main ideas behind this theory; the complete and 
detailed version can be found in [4,1,8]. This treatment is local in the sense that we take 
the region for the time evolution of the problem to be a compact spacelike region S cross 
the positive time (S x R + ). At the time-like boundary of this region, (dS x R + ) -which 
we take to be noncharacteristic- we assume there exist suitable boundary conditions which 
guaranty the uniqueness of the solutions for given initial data. This is so for linear equations 
without constraints, but the non-linear constrained case -which includes general relativity- 
is far from being complete, although results for nonlinear equations and linear systems with 
constraints are available at the present time. 

The key step to prove uniqueness and existence of solutions to hyperbolic systems is to 
establish an a priori inequality which bounds a certain norm of the assumed solution at some 
latter time by a multiple of the norm of the initial data. This inequality is called the energy 
estimate because for physically relevant linear systems the weakest norm of this type that 
can be used to obtain this inequality, and so to assert existence, is just the square root of 
the energy. This a priori energy estimate is: 

Given a solution Uq which for every t G [0, T] is in H m (S), there exists a constant C such 
that given any other solution u k G H m (S) sufficiently close to the first one we have Q: 



where H m (S) is a generalization of the usual Sobolev spaces of order m, which not only 
include space derivatives up to order m, but also time derivatives and crossed time-space 
derivatives up to the same order. In spite of that, the norm at t = is bounded by the usual 
Sobolev norm of order m if one uses the evolution equations to trade all time derivatives 



1 Here we consider first order systems, and so the initial data is just u(0). 
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appearing there for space derivatives. The bound using just the initial data -and not their 
time derivatives- is the one that is important in establishing existence of solutions; we shall 
come to this again when analyzing systems with different time scales. The minimal value 
of m for which the energy estimate is valid for a generic (non-linear) symmetric hyperbolic 
system -a concept we introduce bellow- is the smallest integer larger than n/2 + 2, where 
n = dim S, for this guarantees that we have pointwise bounds on u, a sufficient step in the 
non-linear case to obtain the inequality. 

A sufficient condition for a system of evolution equations to obey the above inequality is 
that it be symmetric at a solution u k , that is, a system that can be written in the following 
form: 

AV) = A\(u k )V a v? + B t (u k ), (2) 

where the matrix is symmetric and positive definite and the vector matrix A a ij(u h ) is 
symmetric. The above matrices and Bi are supposed to be smooth functionals of the vector 
u k . The connection V a is some arbitrary connection on S. 

If we scale some of the components of A a ij with a factor 1/e, as it happens in systems 
with different time scales, where 1/e is the largest speed of the system, much larger than 
the others, then the constant appearing in the energy estimate, equation ([!]), will in general 
appear also scaled with a 1/e factor and in the limit e —>■ as this factor goes to zero one 
losses the estimate. This is not always the case, and as shown in [4,1,8] the dependence of 
C on e is regular if f\ both the following conditions are met: 

1. ) The matrices A a ij appearing in the symmetric hyperbolic system have the following 
form: 

A a ij{u k 1 e) = —A^ij + A^ij(u k , e), (3) 

with Ayij being constant matrices, and A\^ regular in e. The reason for this is simply 
that in this case the singular terms go away on integration by parts. 

2. ) The matrix A°ij, and the vector Bi are regular in e. 

There is another way to loose the estimate, namely if the norm at t = -i.e. on the 
initial data- blows up in the limit (here we refer to the fact that in trading time derivatives 
of the initial data for their space derivatives one again encounters the singular part of A a ij). 
The only way to avoid this singular behavior is to choose very special initial data, that is 
to constrain the initial data by imposing on them some (elliptic) differential equations that 
guarantee the good behavior of the bound on the limit. This process of selecting particular 
initial data, and so the so called slow solutions, is called initialization. To summarize, the 
above theory tells us that if the symmetric hyperbolic system satisfies conditions 1.) and 2.) 



2 We are using Penrose's abstract index notation, by which latin supra-indices denote vectorial 
type entries of a tensor and sub- indices co- vectorial ones. 

3 Actually there exist more general and so more sophisticated conditions; we do not include them 
here because they are not needed for our purposes. 
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above then there are certain solutions coming from initialized data which depend smoothly 
on the parameter e, in the sense that we get the a priori estimate (1). 

In the next section we define appropriate variables for general relativity in such a way 
that the resulting matrix A a ij satisfies condition 1.). Unfortunately condition 2.) is in 
general not satisfied, for the resulting vector Bi contain singular (w.r.t e) terms. These 
terms can not be completely eliminated with gauge conditions, but their undesirable effects 
on the energy estimates actually can. In section 4 we show that all singular terms on the 
integral of UjB i {uk)5 : ' 1 , and similar terms appearing in higher order energy expressions, can 
be annihilated by a particular gauge choice. This gauge fixing is given by elliptic equations 
that must be solved in each time slice, but since these equations depend on the dynamical 
variables we must then treat a coupled hyperbolic-elliptic system of equations. We give 
arguments showing that, nevertheless, the usual estimates can be obtained. 

III. THE SYMMETRIC HYPERBOLIC SYSTEM 

In this section, as a first step to treating Einstein equations as a system with different 
time scales, we cast them as a first order symmetric hyperbolic system. The system we 
present here has this remarkable property for arbitrarily given lapse and shift variables, and 
was found in collaboration with R. Geroch. Further details of that work will be published 
elsewhere. We shall make use of this property when studying the newtonian limit, for this 
requires the choice of a specific gauge, which only arises as a posteriori consequence of 
having the system in a symmetric hyperbolic form. We remark that we were not able to 
obtain a regular limit in the harmonic gauge; the same problem is already present in the 
corresponding limit for the electromagnetic field (in the Lorentz gauge) and has to do with 
the coupling to the fluid. 

A. Variables and Equations 

We consider now the dynamical problem of general relativity, namely, the temporal 
evolution of a 3-dimensional spatial metric. As usually, we take a (3 + 1) decomposition of 
spacetime, so we take it to be foliated by spacelike surfaces S t , which are the level surfaces 
of a function t. The normal n a to these surfaces is Ndt for some function N. The spacetime 
metric g a b induces a metric q a b on the spatial slices: q a b = g a b + n a n b . To this variable 
corresponds a canonically conjugate one, its momentum: n ab = ^{q ab K c c — K ab ), where 
K ab = ^£ n cq ab is the extrinsic curvature of the surfaces and y/q is the square root of the 
determinant of the 3-metric 

In the variables (q ab , 7i ab ) Einstein equations split into a set of evolution equations and a 
set of constraints on the initial data (see for instance [9] []): 



4 The square root of the determinant of the metric is just a shorthand for the volume element of 
that metric 

5 The definition of momentum given in this reference has a difference with the momentum assumed 
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if b = -E- ( 2 fc ab - <f 6 tt c c ) - 2D (a m (4) 

£ 2 n ab = -v^iV {R ab - ]^Rq ab + -^(q ab D c D c N - D a D b N)} 

Nf 2 1 1 

+— {fc ab 7v c c - 2n a c n cb + -q ab (n cd n cd - -vr c c ^)} 

+2S ab e 4 + e 2 N c D c fc ab - 2e 2 7T c(a D c N b) (5) 

R - -(n cd n cd - U c c 4) ~ 4e 2 p = (6) 
q a 

- 2D c 7r ca = 4e 2 J a (7) 

Here, N a is the shift vector; (') = £ t a [t a = ^-n a + N a ) and D a is the derivative operator 
on the slices associated with q ab . The parameter e is the inverse of the speed of light, which 
here will be taken to be the fastest speed when the above system is considered a system with 
different time scales. The way in which this parameter appears on the equations, at this 
level, is determined by dimensional considerations, see the appendix for the rules by which 
we assign dimensions to the different tensor quantities. 

Note that in these equations, second derivatives of the metric are involved, for they 
appear in R ab . 

If equations (6), and (7), the so called constraint equations, are satisfied at any given 
instant of time -i.e. by the initial data- then equations (4) and (5) imply they are satisfied 
at all times. 

The tensor fields p, J a , and S ab that appear in the equations are the different projections 
of the energy- momentum tensor of the matter fields on S t . We are assuming in what follows 
that the matter fields obey by themselves symmetric hyperbolic equations which are coupled 
to the metric in such a way that after we cast Einstein equations in a symmetric hyperbolic 
form, then the whole system would be symmetric hyperbolic. We therefore ignore the matter 
equations. 



B. Conformal Transformation and Lapse Shift Scaling 

Let q a b = N 2 q ab and therefore q ab = N~ 2 q ab . This choice of conformal factor leads us to 
the following expression of the Ricci tensor of the conformal metric: 

Rab = Rab - ^(q ab D c D c N + D a D b N) + ^-D a ND b N (8) 



here, but, having this in mind, the set of the equations we present can be obtained from the 
equations given in this reference by a straightforward calculation. 
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Thus, in the evolution equation for tt second derivatives of the lapse function appear 
only in the form of the laplacian q ab D a Df,. This is highly convenient since, as we shall shortly 
see, this term can be eliminated in favor of the mass density. 

Let 

where U is an arbitrary function. This choice will help us to get rid of some combinations 
of second derivatives of the metric that hamper the system to become symmetric hyperbolic. 
The choice {7 = would correspond to the temporal harmonic gauge, see for instance [5]. 
Here we don't restrict ourselves to this gauge, for U remains arbitrary. 

Nevertheless notice that, up to first order in e, we do have the temporal harmonic gauge. 
The 2 nd order correction to it, U, will be latter identified with the newtonian potential. This 
is in agreement with Nester and Kiinzle [5]. 

Second derivatives of the lapse will then be proportional to second derivatives of the 
newtonian potential plus other terms involving second derivatives of the metric, which will 
be arranged by means of a redefinition of variables in the next section. In this way, there 
will be no second order derivatives of the newtonian potential other than the laplacian. 

The same happens with the curvature scalar [], so the scalar constraint is: 

N 2 R + t f d D c D d N - l-q cd D c ND d N - Ae 2 p - -{^ cd - ~7f^) = (10) 
JM JS Z q I 

So far, the lapse has been redefined and the conformal factor chosen so that a new 
variable U appears conveniently for latter purposes. There are no restrictions yet on this 
variable, which doesn't need to be fixed in order to set a well posed initial value formulation, 
i.e. to get a symmetric hyperbolic system. Restrictions only will come by means of a gauge 
fixing procedure to obtain a regular newtonian limit, thus justifying the association of U 
with the idea of a newtonian potential. 

For similar purposes we also re-scale the shift, 

N a = e 2 N a (11) 



C. New Variables 

To get the symmetric hyperbolic system we define the following variables: 

r ab c = ^- 3 (d c q ab - \q ab q ed d c q ed ) = ^f c {^qq ab ) (12) 



6 We mention this in order to have a complete description, although this form of the scalar con- 
straint is not relevant to the actual computations leading to the final form (equation 16). 
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(13) 



Here, d c is the derivative operator associated to a flat e ab . For the purpose of studying 
the Newtonian limit we introduce here explicitly the inverse of the speed of light (e) on the 
definition of the variables and formulae. This is of no relevance for obtaining the symmetric 
hyperbolic system. 

It is expected that y/qq ab [] will differ from the above flat metric in order e 3 . We shall 
assume, and in fact then assert, that y/qq ab = y/ee ab + e 3 h ab , with h ab e-smoothQ. Therefore, 
derivatives of the metric density will be of order e 3 . One can check that this is so, for 
instance, in Schwarzschild. 

We make use of the constraints to rearrange terms in the evolution equations. With 
appropriate factors, we add the scalar constraint to the equation for it ab to eliminate the R 
term, and we add the vector constraint to the equation for the new variable r ab c . The new 
terms added in this way are necessary to symmetrize the system (p ab , r ab c ). This, of course, 
involves the appearance of extra source terms in the evolution equations: 

p ab = - q iN^(q cd d d r ab c - 2q< a d c r h)d d ) - 2q^q ab (AU - p) 

+2S ab + £ 2 N c dcp ab + £F ab( £ ^ r< fe ^ p de ? qjj^ (14) 



. 1 



f ab r 



N ~(d c p ab -25i a d d p^ d ) 



q 



i e 



+ \q ab d c d d N d - q d ( b d c d d N a ^} 

+ ^-5i a J b) + -q ab d c {NU) + e 2 N d d d r ab c 
q<L e e 

-2eNq ab N d d c d d U 

+eF ab c (e, r ab c , p ab , d c U, U, d c N a ). (15) 

Here A = q ab d a d b , F ab = F ab (e,r ab CJ p ab ,d c U,d c N a ), 
F ab c = F ab c (e, r ab c ,p ab , d c U, U, d c N a ), and all other F's appearing in the equations from now 
on are smooth pointwise functions of all their arguments. 

Of course for the solutions to this system to be solutions of Einstein equations one needs 
to impose on the initial data the usual constraints, which in these variables take the form: 

AU - \ed c r cd d = p + e 2 F(e, r ab C} p ab , d c U) (16) 
- 2d c p ca = 4J a + e 2 F a (e,r ab c ,p ab , d c U), (17) 



7 Notice that this metric density agrees with the dynamical variables of M. Lottermoser in (7). 
this will follow directly once we show that r ab c is pointwise bounded. 
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and the extra one: 



^=^^(7^), (18) 

which if satisfied initially then, -modulo boundary conditions-, hold at all times. To see 
this take time derivatives of the constraints and using the above equations get a symmetric 
hyperbolic system. Thus, if the boundary conditions for the original system, eqn's (4, 14, 15, 
16, 17, 18), are such as to ensure uniqueness of solutions for this derived system, then the 
constraints would remain zero for all times and the original system is equivalent to Einstein's 
equations, (eqn's (4,5,6,7)). It is interesting, and perhaps disturbing, that the constraints 
propagate at a different speed than (but numerically proportional to) light. 

We claim that the above system is symmetric hyperbolic, for any given choice of lapse 
and shift (U, N a ). To see this one can compute the resulting matrices A°ij, A a ij, and check 
their symmetry. We faind it convenient to split the matrices into pieces using compounded 
subindices. Thus we define, for instance, A^ r that part of the matrix A^ that acts on p ab 
and has image in the space of the r a6 c 's. We take as the matrix A^, one having only the 
following nonzero components: A° qq = S^S b , (we reserve the first letters of the alphabet for 
contraction with vectors in the domain and the last ones for the image), A^ p = q~ x S^Sl, 
A° rr = S^5 b 5^.. The only nonzero components of the resulting A\ - vector valued matrix are 
then: 

Kr = -^f-m^l - 2S c (a 8l/ q ), and 

From these formulae by contraction with two arbitrary vectors, it is easy to see the 
symmetry of A\j. Alternatively, one can compute the time derivative of the first norm 
appearing in the energy estimate ( |1|), i.e. the energy, and see that after integration by parts 
it only depends on the dynamical fields (y/qq ab , p ab , r ab c ), and not on their derivatives. The 
expression for this energy is: 

E(t) = \l {r ab c r ab c + q~ W} dZ, (19) 

where we have raised and lowered indices using the conformal metric, q ab . 
Its time derivative is: 

E(t) — — J \r ab c r ab c + q^p^Pab] dY, + algebraic terms in the fields. (20) 

Using equations (|T4] ) (pT5|), and integrating by parts, we obtain: 

E{t) = - [ \-r ab c {q ab d c d d N d - q d{b d c d d N a ^) 

4ivi 1,2 

-r ab a J b - -r ab c q ab d c (NU) 



q4 £ £ 
-2eNr ab c q ab N d d c d d U 

-2q-^p ab q ab (AU - p)\ dY 

+ e regular algebraic terms in the fields and in 
first derivatives of U and N a . (21) 
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We have explicitly written the terms which contain second derivatives of the shift and 
of the potential, or are singular in e, for latter purposes, although they have no relevance in 
this section because here we are still assuming the lapse and shift are freely given, and we 
are not yet taking the limit e — > 0. 

This symmetric hyperbolic system, with sources having symmetric hyperbolic evolution 
equations (such as a perfect fluid), is sufficient to assert the existence of solutions to the 
relativistic equations for any given lapse and shift, -i.e. for any gauge-, and any given value 
of the parameter e, as small as desired but different from zero. In what follows we are 
not going to be interested in existence, since it can be derived from this system or others 
systems in the literature. We are going to be mainly concerned with the e-smoothness of 
these solutions for which we assume existence. 



IV. THE NEWTONIAN LIMIT 

The singular behavior in e of the time derivative of the energy prevents us from obtaining 
an energy estimate, equation (|TJ), with constant independent of e, and therefore, as discussed 
in the second section, we do not control the behavior of the solutions in the limit e — > 0. 
From the viewpoint of the general theory described on the second section the problem arises 
because in our case the Bi{u k ) on equation (2) are in fact singular with respect to e. These 
singular terms can not be eliminated from equation (2) by any gauge condition, but as 
we shall see the singular terms that they generate on the energy estimates can indeed be 
eliminated by choosing a convenient gauge. 

To see this we do further integration by parts on the expression for the time derivative 
of the energy to obtain: 

1 ANi 

E(t) = - I { -r ab a { — —J b - q cd d c d d N b ] 



s t I e q~ 

-- £ r ab c q ab d c {2NU + d d N d ) 

-2sNr ab c q ab N c AU - 2q-^p ab q ab (AU - p)| rfS 

+ e regular algebraic terms in the fields and in 
first derivatives of U and N a . (22) 

We now use the scalar constraint equation to get rid of the laplacian of U, and get: 



E{t) = - M - r °*> c 



s t e 



^Lljb _ q ^d c d d N b - e 3 Nd b r de c q de N c 



Q' 2 

-q-^d b p de q de ] ~r ab c q ab d c (2NU + d d N d )\ dZ 

+ e regular algebraic terms in the fields and in 
first derivatives of U and N a . (23) 
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We observe that there are two different types of singular terms, one is a factor of r a 
the other of r ab c q a b. To eliminate the first term we choose a gauge (i.e. a selection of lapse 
and shift) that makes r ab a = 0. But to achieve this we need U to satisfy: 

AU = p + e 2 F(e,r ab c ,p a \d c U), (24) 

since otherwise the constraint equation (16) would imply that dbr ab a is different form 
zero. We then adjust iV a such that r ab f ) = 0. This implies the following equation for N a : 

d c D c N b - d b D c N c - 4{J b + d b NU) + 2e 2 Nq cb N d d c d d U - e 2 G b = 0. (25) 

A necessary condition for this equation to have solutions is that the above equation for 
U be satisfied. In fact, from Bianchi identities it follows that the divergence of the above 
equation vanishes identically, regardless of N a . Equation (25) does not fix N a completely; 
in fact as it stands it is not even elliptic. We use the remaining freedom in the lapse to 
choose 

D d N d = -2NU, (26) 

getting rid of the second term. To see that this is possible we add to equation (25) the 
term q bc d c (4NU + 2D d N d ), thus eliminating the terms with U from it and also rendering it 
elliptic: 

d c D c N b + d b D c N c - AJ b + 2e 2 Nq cb N d d c d d U = e 2 G b (e, r ab c ,p ab , d c U). (27) 

Since the divergence of equation (25) vanishes identically the divergence of equation (27) 
must then be: 

= A(ANU + 2D d N d ). (28) 

Uniqueness of solutions to Laplace's equation then implies equation (26). 

We have thus eliminated from the time derivative of the energy all e-singular terms. 
The price we pay for this regularization is that now we must consider a mixed symmetric- 
hyperbolic-elliptic system of equations, the hyperbolic part being equations (4, 14, 15), the 
elliptic part being equations (24,27) []. The initial constraint Q equations are now (17,18) and 
r ab a = 0. The time derivative of U can be eliminated from all the equations by using equation 
(26). The authors do not know of any general treatment of mixed symmetric-hyperbolic- 
elliptic systems, and so in what follows we present an argument leading to establish an 
energy estimate, and thus existence of solutions for mixed systems. Because the scope of 



9 Equations (24,27) form an elliptic system for (U, N a ) as long as the tensor q ab — e 4: N a N b is 
positive definite. This in turn is a condition (for e ^ 0) on how big J a is. 

10 We call constraint equations those which only involve initial data, and which once they are 
solved at the initial slice they remain valid for all later times by virtue of the evolution equations. 
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this argument is broader than the topic here treated, and it involves further techniques, we 
shall give the details elsewhere. 

The argument is as follows: If we had Garding estimate for the elliptic variables in 
term of the hyperbolic ones, following the procedure in [1] with only slight modifications 
we would arrive to an energy estimate for just the hyperbolic variables -in a sense the 
elliptic variables are just non-local (but smooth) functionals of the hyperbolic ones-, and 
the existence of solutions would follow in the usual way Thus we just need Garding's 
estimate for the elliptic fields, and for that we need to make sure the elliptic equations are 
injective at each time slice. We consider, for given initial data for the hyperbolic fields, the 
elliptic system as depending on two parameter families, e, and t, the evolution time. It is 
easy to see that for e = the elliptic system -with suitable boundary conditions- has a 
unique solution, regardless of what the hyperbolic fields are. Since the space of injective 
operators is open, then given any initial data there will exist an Eq such that for all e 6 [0, Eq] 
the elliptic part of the system is injective at t — 0. But for the same property, and since we 
are assuming smoothness on t (through the a priori smooth dependence of the hyperbolic 
fields on t), there will exist a To > such that the elliptic part of the system will remain 
injective for all times smaller than or equal to To. Thus we do have Garding's estimate 
during a finite period of evolution. 

We have established an a priori energy estimate for general relativity which remains 
valid even at the e = limit. We remark that the matrix A^- that appears in our equations 
is regular in e times those dynamical variables whose equations are e-singular, and so as 
remarked in [1], it is possible to get estimates which only involve spatial derivatives of the 
fields in the energy norms PI Therefore the initialization process that we carry out bellow 
is not needed as far as to get a priori estimates for norms involving only spatial derivatives. 
Nevertheless, some initialization is needed to obtain regularity in time of solutions, for in 
the singular equation case one can not use the equations to conclude smoothness in time 
from smoothness in space. 

V. INITIALIZATION 

For the reason given above, it is important to control the time derivative of the dynamical 
fields, that is, to get a priori estimates for energy expressions which include L 2 norms of 
these time derivatives. This gives boundedness in time of the family of solutions. 

We now study the conditions on the initial data for these energy norms to be initially 
bounded; that is all what is needed as made clear by the results obtained in the last section. 
Obviously the terms having no time derivatives of p ab and r ab c do not have any singular 
dependence on e. We first encounter problems with the terms p ab p a b + r ab c i r ab c , and therefore 
we must choose the initial data -to initialize it- so that they are regular. Writing the initial 
data as, U\ = U + eU\ + etc. we have that the condition for r ab c to be regular at the 
initial time is: 



That is, in this case it is possible, by using the evolution equations, to trade time derivatives by 
spatial derivatives of the variables without getting e-singular contributions. 
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d c p<*> = -e d ^ h d c d d N*\ (29) 

that is, 

pt h = -e d( > h d d N$. (30) 

Notice that this condition is consistent, since using the zeroth order part of equation 
(27), we get: 

daPo = -\e db d d d a N« - l -e da d d d a N h = -2j\ (31) 

that is the zeroth order term of the vector constraint, equation (17). 
The condition for p ab to be regular is that 

e cd d d rf c = 0, (32) 

which implies that initially we must take Wq = 0. 

If further smoothness in time is required for some application, then one has to consider 
further time derivatives in the energy norm, and continue the initialization procedure. 

In the compact hypersurface case further orders give rise to equations that can be solved 
at each order without apparent obstructions, so in that case one would conjecture infinite 
smoothness in time. Existence of solutions to this finite hierarchy of equations and of the 
constraints for the asymptotically flat case deserves further study. See references [6,7] for 
different approaches to get these conditions, and the existence of solutions to them. 



VI. CONCLUSION 

We have shown the existence of an a priori energy estimate which holds for any value 
of the parameter e, including the zero. These estimates are an important step towards the 
proofs of existence^] and smoothness of slow solutions, i.e. solutions satisfying the initial- 
ization conditions, for any given matter source system which is itself symmetric hyperbolic, 
which is regular in the gravitational variables we have used, and which has a non-relativistic 
current. In particular this is the case for perfect fluid sources with any non-relativistic initial 
configuration. 

Which setting can we embed our calculation in to obtain complete results on the exis- 
tence and e-smoothness of near Newtonian solution to Einstein's equations? It is clear that 
the main obstacle for that is to find a setting where terms either make a negative or null 



12 Notice that our final (regular) system is a mixed one: symmetric hyperbolic for the dynamical 
fields variables, elliptic for the lapse and shift fields. In this case to establish the a priori energy 
estimate one needs Garding's a priori estimate for the elliptic part. 

13 To be precise, besides the a priori estimate to establish existence one needs an approximating 
sequence of trial functions. 
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contribution to the energy estimate, as to justify the fact that we have ignored them. The 
simplest case where this is so is when there is no boundary, as would be a spacelike closed 
cosmology admitting a Newtonian cosmology limit. This would also be the case if we had 
an initial-boundary value formulation of general relativity admitting maximally dissipative 
boundary conditions. Unfortunately this problem, which is badly needed for making confi- 
dent numerical computations, is not solved yet, partially due to the presence of constraints, 
and partially due to the fact that very little is known about the initial-boundary problem 
for nonlinear hyperbolic systems.^ For this case -we conjecture- one would find uniformly 
locally smooth one parameter families of solutions for any given (small enough) Newtonian 
solution. Another way to state this conjecture would be the following: Given any Newtonian 
solution and any integer m, there exists a general relativistic solution which stays near 
the Newtonian solution to order e m for a Unite time interval. The thir case of interest is the 
asymptotically flat one where again we do not have to worry about boundary conditions. 
The main difficulty here is that one can not use the impressive machinery of Weighted 
Sobolev Spaces and their corresponding Sobolev inequality for unbounded domains, for the 
use of radial functions as weights would introduce unavoidable ^-singular terms on the cor- 
responding energy estimate. Thus, since the Sobolev inequality is badly needed to handle 
non-linearities, we have to resort to Sobolev spaces which also include some given number 
of time derivatives and use the boudedness in time of the evolution region to get the desired 
inequality. But this implies that we need initially some smoothness in the time direction for 
our Sobolev norms to start finite, that is we need a grater degree of initialization. Thus, it 
seems to be the case that e-smoothness in the asymptotically flat case is tied to the absence, 
or presence, of obstructions to solve the hierarchy of initialization equations to the needed 
order. 

One could raise the question of why we need to impose a gauge condition to get this 
results, since after all the theory is gauge invariant. A partial answer to this is that to 
establish the a priori estimates we have obtained one has to treat the fields in some fixed 
gauge, since it is easy to make diffeomorphisms which are singular with respect to e in the 
limit e —>■ 0. But it seems to us that still our gauge conditions could be relaxed somehow, 
for instance it should be enough to fix them up to some order in e and not to all orders as we 
have done We believe that the only gauge (up to first order in e) in which the equations 
are regular is the above one, although we do not have a proof of that. Of course one can 
pretend less, that is smoothness of ep ab , and of er ab c , in that case the only non-singular term 
in equation (22) is the one proportional to (AU — p), and so the gauge: AU = p, N a = 
suffices, but this is not optimal. 



One of the authors has already some results for the linearized Einstein's equations in the har- 
monic gauge. 

15 Near in the sense of the energy norm. 

16 Note that we need to fix it to all orders to ensure that our system is symmetric hyperbolic- 
elliptic. 
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There are other issues which deserve further study. 

For analytical and numerical studies of the characteristic problem (the initial time for- 
mulation along a future light-cone), it is of interest to treat the related Newtonian limit, see 
for instance [10], and establish similar results to the ones here obtained. In particular this 
should be important as another justification to pick initialized data as data with little extra 
radiation apart from the one coming from the matter sources. For related results see [3] . 

Once one has control of the incoming radiation not generated by the sources one can 
start to consider in a rigorous way the back reaction of radiation on the sources. There 
are methods available to treat those effects along the lines we have considered here, and we 
think they deserve some attention. 

APPENDIX A: DIMENSIONS 

In this appendix we briefly introduce the concept of dimension or of units for geometrical 
objects, and give the rules we followed for assigning dimensions to several geometrical objects 
that we used in the third section of the paper in order to be able to determine where e 
appears on the equations. Let M be a manifold, and C(M,R) the algebra of smooth real 
functions on it. We introduce dimensions by enlarging this algebra to the cartesian product 
of C(M,R) with a discrete abelian group, G. Each element of this new algebra is then a 
pair (/, T), where / G C(M, R) and T G G is the dimension of /. The product is the usual 
one: (/, T) x (g, L) = (fg, TL) and the sum is only defined within pairs with the same group 
element: (/, T) + (g, T) = (/ + g, T). We shall, as usual, omit the pair and only write / for 
(/, T), and we will denote the projection to the group entry as [(/, T)] = [f] = T. 

We shall require the coordinate functions to have a dimension different than unity, say L, 
and define the dimension of a vector n a to be the element [n a ] of G such that when the vector 
acts on any function, / G C(M,R) x G, it gives a function of dimension [n(f)] = [n°]^[/]. 
This definition is equivalent to assigning the same dimension to the components of the vector. 
Note also that the Lie bracket of two vector fields yields a vector of dimension equal to the 
product of their respective dimensions divided by the coordinate functions dimension. 

We define the dimension of a covector m a to be the element [m a ] of G such that when 
the covector acts on any vector, n a it gives a function of dimension [m a n a ] . We extend these 
definitions to tensor fields in the obvious way. Note that then the connection has dimension 
of L _1 and the Riemann tensor has dimension of L~ 2 . If a metric tensor is present, then 
to be consistent with the formulas in a coordinate system, the metric -and therefore its 
components- has to have as dimension the group identity, [g a b] = E. With this convention 
then the length of a vector has the same dimension as the vector, and raising and lowering 
indices do not change the dimension of the objects. 

For our application we assign: to the function t that foliates the space-time the dimen- 
sion [t] = T, to the unit normal to the surface, [n a ] = E, to the time flow vector [t a ] — ^, (so 
that [t a V a t] = E), to the inverse of the speed of light, [e] — j j . The dimension of any other 
object is defined following the rules stated above, in particular; [q ab = g ab — n a n b ] = E, since 

t a = | n a + i ya ? ^ = ^ &nd [^a] = L ? = 2 V (a n b) ] = \, [lT ab = &(q°» R C c ~ K ab )\ = ± 
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